Dynamic criticality in glass-forming liquids 
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We propose that the dynamics of supercooled liquids and the formation of glasses can be under- 
stood from the existence of a zero temperature dynamical critical point. To support our proposal, 
we derive a dynamic field theory for a generic kinetically constrained model, which we expect to 
describe the dynamics of a supercooled liquid. We study this field theory using the renormalization 
group (RG). Its long time behaviour is dominated by a zero temperature critical point, which for 
d > 2 belongs to the directed percolation universality class. Molecular dynamics simulations seem 
to confirm the existence of dynamic scaling behaviour consistent with the RG predictions. 
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The recipe for making a glass is simple [l| : rapidly cool 
a liquid through its melting point to avoid crystalliza- 
tion. Cool it further, and the liquid eventually becomes 
so viscous that it forms a non-crystalline solid, or glass. 
Glasses are common in nature, but our theoretical under- 
standing of their formation is poor . Here we offer an- 
alytical and numerical results to support the proposition 
that the dynamics of glass-forming supercooled liquids is 
controlled by a zero temperature dynamic critical point. 

Our starting point is the real-space description of su- 
percooled liquids studied in 0, PI and based on ideas and 
models originally proposed in jjj, IE • In this approach 
the phenomenon of dynamic heterogeneity [HIEE! plays 
a central role. The link between dynamic heterogeneity 
and glass formation is the subject of much current re- 
search. If this link is verified it will be an indication 
that the slow dynamics of glass formers is governed by 
dynamic spatial fluctuations (see [llj, ll2j for alternatives 
based on thermodynamics), in contrast with the assump- 
tion of homogeneity of mode coupling theories 0, Q] . 

The microscopic coarse-grained approach of 0, 0] re- 
lies on two observations: (i) at low temperature very few 
particles are mobile, and these mobility excitations are 
localized in space; (ii) mobile regions are needed to allow 
neighbouring regions to themselves become mobile. This 
is the concept of dynamic facilitation [5|,l3- We show that 
this picture can be cast as a dynamical field theory and its 
scaling behaviour derived from a dynamic renormaliza- 
tion group (RG) analysis. We find that scaling properties 
are determined by a zero temperature critical point. For 
the simple case we consider, and for d > 2, this critical 
point is that of directed percolation (DP) [l^J- We also 
show, by performing extensive molecular dynamics simu- 
lations, that supercooled Lennard- Jones binary mixtures 
display scaling behaviour consistent with that predicted 
by RG. 

We build an effective microscopic model as follows. A 
supercooled fluid in d spatial dimensions is coarse-grained 
into cells of linear size of the order of the static correlation 



length as given by the pair correlation function. Cells are 
labeled by a scalar mobility field, rii, identified by coarse- 
graining the system on a microscopic time scale. Mobile 
regions carry a free energy cost, and when mobility is low 
interactions between cells are not important. Adopting 
a thermal language, we expect static equilibrium to be 
determined by the non-interacting Hamiltonian 0, 

N 

H = J2"i- (!) 

i=l 

At low mobility, the distinction between single and mul- 
tiple occupancy is probably irrelevant. We assume the 
latter case for technical simplicity. The dynamics of the 
mobility field is given by a master equation, 

d t P ({n} ,t)=J2d (W) A P ({n} , t) , (2) 

i 

where P ({n} ,t) is the probability that the system has 
configuration {n} at time t. The local operators Ci en- 
code the existence of local quanta of mobility. For non- 
conserved dynamics they describe the creation and de- 
struction of mobility at site i, 

CiP{m,t) = 7(n i + l)P(n i + l ) t) + pP(n i -l,t) 
-( 1 n i + p)P(n i ,t), (3) 

where the dependence of P on cells other than i has been 
suppressed. The rates for mobility destruction, 7, and 
creation, p, are chosen so that @ obeys detailed balance 
with respect to at low temperature, p/7 = e _1 / T w c, 
with c = {rii); the brackets indicate a thermal average. 
The average concentration of excitations c is the con- 
trol parameter of the problem. Ci ({n}) is the kinetic 
constraint designed to suppress the dynamics of cell i if 
surrounded by immobile regions. It cannot depend on rii 
itself if © is to satisfy detailed balance. To reflect the 
local nature of dynamic facilitation we allow Ci to depend 
only on the nearest neighbours of i and require that 
Ci is small when local mobility is scarce. 
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The large time and length scale behaviour of the model 
defined by Eqs. CJ-© is derived from the analysis of 
the corresponding field theory. The technique to re- 
cast the master equation (0 as a field theory is stan- 
dard 0, 0], so we only outline the procedure. One 
introduces a set of bosonic creation and annihilation op- 
erators for each site i, a\ and at, satisfying [a\, cij] = Sij 7 
and passes to a Fock space defined by a state vec- 
tor |*(t)) ee E{n«}P('*i. n a. • • • ,*)4™4 n2 ...|0>. The 
master equation then assumes the form of a Eu- 
clidean Schrodinger equation, d\^(t))/dt = — H\^f(t)), 
with H = J2i^i({ a ] a j})^i°^ ■ The unconstrained piece 
H- ^ reads 

Hf ] =- 1 {a l -a\a l )-p{a\-l). (4) 

The Hamiltonian (0J describes the creation and destruc- 
tion of bosonic excitations with rates p and 7. The evolu- 
tion operator e~ Ht can then be represented as a coherent 
state path integral [l7] weighted by the dynamical action 

«S[#,<M„] =J2j ° dt [tidtfo + Hitt*,®] , (5) 

where we have suppressed boundary terms coming from 
the system's initial state vector. The fields 4>*(t) and 
<pi(t) are the complex surrogates of a\ and a,, respec- 
tively, while Hi has the same functional form as (0J with 
the bosonic operators replaced by the complex fields. At 
the level of the first moment we have (rij) = (<fo), so we 
may regard 4>% as the complex mobility field. The last 
step in the passage to a field theory is to take the con- 
tinuum limit, J2i ~^ a~ d f d d x, (f>i{t) — > a d (j)(x,t), and 
4>*(t) — * (f)*(x,t), where a is the lattice parameter. 

The definition of the model is completed by specifying 
the functional form of the kinetic constraint. The sim- 
plest non-trivial form is the isotropic facilitation function, 
C{ = ^rij, where the sum is over nearest neighbours of 
site i. With this choice we expect the one-spin facilitated 
Fredrickson- Andersen (FA) model in d dimensions 
to be in the same universality class as our model. Dif- 
ferent choices for the operators and C lead to field 
theoretical versions of more complicated facilitated mod- 
els. A diffusive H^°\ for example, would correspond to 
a constrained lattice gas like that of Kob and Andersen 
0,0|; an asymmetric C, to the East model 0,0 and its 
generalizations. The dynamic action finally reads 

S$, 0, t ] = J d d x J*° dt [$(d t - D Q \7 2 - 4 m) ) <jy 

+#(A 1) + i/« V 2 )0 + <MA 2) + ^ 2) V 2 )# 
-^(z^+ctoV 2 )^) , (6) 

where we have omitted higher-order gradient terms, and 
suppressed boundary contributions coming from initial 



and projection states. We also made the shift (f>* = 1 + 
0, and defined A 1,2) ee 2da d 1 , k^ >v) ee 2dp, v^' 2) = 
ja d+2 , Co = a 2 p. We defined Do = <jq to emphasise the 
emergence of a diffusive term, although in the unshifted 
model there is no purely diffusive process. The action © 
is the starting point for our RG anal ysis . We will leave 
the technical details to a later paper jl9( and here state 
only the most important conclusions. 

The action 10 has the form of a single species branch- 
ing and coalescing diffusion-limited reaction with addi- 
tional momentum dependent terms |l5|. By integrating 
out the response field, we obtain a Langevin equation for 
the evolution of <fi. This equation has a critical point at 
c = 0, i.e. T = 0, describing the crossover from an expo- 
nential decay of mobility at finite c, i.e. T > 0, to an alge- 
braic decay at c = 0. In the absence of noise (correspond- 
ing to neglecting terms quadratic in <p) we find that ® 
admits the Gaussian exponents (z/[, z^, /3 G ) = (5,1,1). 
Here v± and v\\ control the growth of spatial (£j_) and 
temporal (£||) length scales near criticality, £j_ ~ c~ u± 
and £|| ~ c~ v w , while j3 governs the long-time scaling of 
the density, n ~ '. 

These Gaussian power laws are modified by fluctua- 
tions which we treat using the RG. Identifying the un- 
physical microscopic singularities present in © as a con- 
sequence of taking the continuum limit, one invokes scale 
invariance and dimensional analysis to extract the macro- 
scopic scaling; see Fig. ^ We find the following. 

(i) The critical point remains at c — 0. There is thus 
no finite temperature phase transition. For finite T, and 
at asymptotically long length and time scales, the system 
will therefore exhibit Gaussian power laws. 

(ii) Dimensional analysis shows that the upper critical 
dimension of the model is d c = 4. For d < 4 we ac- 
count for fluctuations by studying the behaviour of the 
effective couplings. These couplings come from diagrams 
like in Fig. ^p. For 2 < d < 4 we find that only the 
dimensionless coupling x = AaK^"' X^ p~ e / D 2 matters, 
where A d = 2 2 - d ir- d / 2 T{Z - d/2), e ee 4 - d, and p 
is an external momentum scale. The system exhibits 
scale invariance at the fixed point x* = 2e/3, and all 
other effective couplings are irrelevant, both at this and 
at the Gaussian fixed point, so that all interaction terms 
in © except for X^ and k v are irrelevant. It follows 
that our system is described for intermediate length and 
time scales by the directed percolation (DP) critical point 
and its associated power laws 0]. To order 0(e) they 
are (i/±,v\\,0) = (§ + ^, 1 + ^, 1 - f ). 

(iii) For d < 2 one cannot access the system's be- 
haviour by way of a perturbation expansion around d = 4 
because newly relevant terms in the action introduce di- 
vergences to all orders in perturbation theory. We do not 
expect DP behaviour, therefore. 

Physically, RG identifies the relevant microscopic in- 
teractions in d = 3 and shows that the action © re- 
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FIG. 1: Elements of the RG calculation, (a) Gaussian prop- 
agator and vertices in the frequency-wavevector domain cor- 
responding to the terms in the dynamic action JHJl; time runs 
from right to left, (b) Structure of diagrams dictating the 
effective couplings, (c) A diagram contributing to the renor- 
malization of the effective coupling x. 



sembles that of DP in the case where the coupling con- 
stant for particle self-destruction vanishes. This is why 
the critical point is at T = 0, so that finite T dynam- 
ics takes place in the active phase of DP ^]. More- 
over, simulations of the 3d FA model also confirm the DP 
picture ^!|. For example, for the correlation timescale 
we expect (Dt) ~ (cr) ~ c~"^ , so that t ~ c~ 5 with 
8 = 1 + 2/|i = 2 + e/12. The numerics indeed shows 
that 6 « 2.10 UM, rather than the naive estimate S = 
1 + 2/dw 1.66 U 

Our analysis implies that the slowdown is a dynam- 
ical critical slowing down as the critical point is ap- 
proached from above. Correlation time and length scales 
grow as inverse powers of c. Thermal activation results 
from c ~ e~ x / T . Dynamical scaling is predicted to oc- 
cur when the dynamic correlation length becomes ap- 
preciably larger than the lattice spacing. This happens 
therefore for temperatures lower than TL the onset tem- 
perature for dynamic heterogeneity jj, l2l| . 

Our analytical results apply to systems with isotropic 
facilitation. We therefore expect them to apply to strong 
liquids, and to those which exhibit a crossover from frag- 
ile to strong behaviour 0. While we do not expect 
DP behaviour for the case of anisotropic constraints or 
conserved order parameters [l^, we do expect a zero- 



temperature critical point to be the gene ric feature of 
both strong and fragile glass formers |l9l Eof . We ex- 
pect that the scaling properties of liquids in their fragile 
regime can be described by the field theory for a facili- 
tated model with directional persistence, such as the East 
model 0. 

The field theory suggests the following physical pic- 
ture. The viscosity of a supercooled liquid increases 
rapidly as T is lowered, because the dynamics becomes 
increasingly spatially correlated. A glass is obtained 
when the liquid's relaxation time exceeds the experi- 
mental time scale. The scaling properties of time and 
length scales and therefore the physical properties of su- 
percooled liquids are governed by a zero temperature dy- 
namic critical point. 

As shown in 0, 0] , this picture accounts quantitatively 



for several observations concerning thermodynamic, dy- 
namic, topographic and spatial properties of supercooled 
liquids. However, the existence of critical fluctuations 
and dynamic scaling remains to be proven [22I E3 | . We 
now present numerical evidence to this end. 

We have performed molecular dynamics simulations of 
a well-characterized model for supercooled liquids, the 
canonical 80:20 binary Lennard- Jones mixture defined in 
Ref. 13 ; further numerical details are found in . Ac- 
cording to the above analysis, power law spatial corre- 
lations develop in the dynamics of a supercooled liquid 
when T is lowered. With this in mind, we have measured 
«Sk(q,T), defined as the spatial Fourier transform of the 
following two-point, two-time correlation function, 



C k (r) = 



(fk(0,r)F k (r,r))-(F k (r,r)) 2 
(F k (r,r) 2 )-(F k (r,r)} 2 



(7) 



where F k (r,i) = £,^(0) - r) cos (k • [r,-(t) - r^O)]) 
is a natural local indicator of the dynamics of the liquid. 
In these expressions rj (t) is the position of particle j at 
time t, and r = r(k, T) is defined in a standard way from 
the time decay of (F k (r,i)). Similar functions replacing 
density correlations by particle overlap or velocity cor- 
relations have been discussed 0, 0, • As predicted 
theoretically, we find that the q dependence of S* k is well 
described by the following scaling form, 



(8) 



for all k, although a precise determination of 77 has proven 
impossible because it is small, as is found in the RG ap- 
proach. We are aware of no alternative theoretical pre- 
diction for the correlator JJJ. The scaling function S(x) 
behaves as S(x — > 0) = const, S(x — > 00) ~ x v ~ 2 . More- 
over, it is universal in the sense that it appears to be 
independent of k for a range of inverse wavevectors from 
the interparticle distance up to the size of the simulation 
box. From JSJ, we conclude that universal power law cor- 
relations develop in supercooled liquids. Imposing n = 0, 
we can use Eq. (|8J) to numerically extract the length scale 
£ k (T), which grows when T decreases, as in [2^. 

The above RG analysis predicts scaling laws for length 
and time scales expressed as a function of c, the con- 
centration of mobile regions. This quantity is difficult 
to measure in simulations. However, one can eliminate 
c in favour of I ~ t 1 / 2 , which defines a dynamic expo- 
nent z = S/v±. Similarly, it is convenient to measure 
% k (T) = Sk(q = 0,T), which is predicted to scale as 
X ~ t x ^ a , A" 1 = (2 — rj)v_\_/5. These two scaling pre- 
dictions are tested in Fig. |2 We find that the power law 
scalings anticipated theoretically are well satisfied. We 
find also that the exponents z and A that we measure 
from Fig.0do not depend on the chosen physical observ- 
able. We note that the numerical values of the exponents 
are in reasonable agreement with the RG predictions, de- 
spite the fact that the simulated liquid is thought to be 
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FIG. 2: Dynamic scaling of the susceptibility (top) and the 
correlation length (bottom) measured for various wavevec- 
tors in the binary Lennard- Jones mixture; (Ar) 2 indicates 
the k —> limit. For each k, we independently measure Xk, 
lit. and Tk for T £ [0.45; 2.0]. The parameters xo and £o are 
determined as the numerical prefactors in power law fits, so 
that the universal dynamic scaling observed in these plots is 
obtained without free parameters. Errorbars reflect observ- 
able to observable fluctuations of the exponents, and do not 
take into account systematic errors. The few points outside 
the power law in the bottom figure are all for T = 2.0, much 
above the onset of slow dynamics, T ~ 1.0. 



fragile. A more precise comparison between theory and 
simulation would require a better finite size scaling anal- 
ysis [2^. Note finally that Q| indirectly predicts a sim- 
ilar scaling between \ an d t although neither the value 
of the exponent nor its relation to spatial structures were 
derived. 

In conclusion, this paper presents further jjj theo- 
retical evidence that the slow dynamics in supercooled 
liquids is governed by a zero temperature dynamic crit- 
ical point at which time and length scales diverge. We 
propose that this critical point is responsible for the ex- 
istence of the glass state. 
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